package GSEA::Pvalue::Permutor ;

use strict		; 
use warnings	;

use GSEA::ES::EScalculator ;

sub constrPermutor($$$$$) #constructeur
{
	my ($classe , $S_size , $L_table  , $Npermutations , $Pparameter ) = @_;
	
	
	if( ! defined($S_size) 			){ die("argument 'S_size' is not set , error");}
	if( ! defined($L_table) 		){ die("argument 'L_table' is not set , error");}
	if( ! defined($Npermutations)	){ die("argument 'Npermutations' is not set , error");}
	if( ! defined($Pparameter) 		){ die("argument 'Pparameter' is not set , error");}
	
	my $this = {
					"S_size"		=> $S_size	,
					"L_table"		=> $L_table	,
					"Nperm"			=> $Npermutations,
					"Pparameter"	=> $Pparameter,
					"permutations"	=> 0		,
					"scores"		=> 0  
				};
	
	bless ($this,$classe); 
		
	return $this;
}


sub _generateRandomDistribution( $ )
{
	(my $this) = @_ ;
	
	my $S_size			= $this->{"S_size"}	;
	my $pt_Ltable		= $this->{"L_table"};
	my @L_table 		= @{ $pt_Ltable } 	;
	my $Npermutations	= $this->{"Nperm"}	; 
 
	if( $S_size > @L_table ){die("can not generate a geneset of ".($S_size)." genes , with a list of ".($#L_table +1)." available genes !");}
	
	my @permutations ; 
	
	for(my $i=0 ; $i<$Npermutations ; $i++)
	{
		# Shuffleing
		for( my $ch = 0 ; $ch< @L_table ; $ch++)
		{
			my $a = int( rand(@L_table) ) ;
			my $b = int( rand(@L_table) ) ;
			
			my $temp		= $L_table[$a];
			$L_table[$a]	= $L_table[$b] ;
			$L_table[$b]	= $temp ;
		}
		
		# creating a geneset , with the first $S_size genes
		my @perm   ;
		for(my $j=0 ; $j< $S_size ; $j++)
		{
			push( @perm ,  @{ $L_table[$j]  }[0] );
		}
		
		push(@permutations , \@perm);
	}
	 
	 $this->{"permutations"} = \  @permutations ;
}

sub _calculateGeneSetScores( $ )  
{
	( my $this ) = @_ ;
	
	my $pt_Ltable	= $this->{"L_table"};
	my @L_table 	= @{ $pt_Ltable } ;

 	my $pt_geneSets = $this->{"permutations"};
 	my @geneSets	= @{ $pt_geneSets };
	
	
	my @scores;
	
	my $counter = 0 ;
	
	# calculating scores :
	print("\n");
	for( my $i = 0 ; $i< @geneSets ; $i++ )
	{
		my @S_table = @{ $geneSets[$i] } ; 		
		
		my $calculator = GSEA::ES::EScalculator->constrEScalculator( \@L_table ,  \ @S_table , $this->{"Pparameter"} ) ;
		$calculator->searchCommonGenes();
		$calculator->CalculateEScourb();
		(my $ESmin , my $ESmax ) = $calculator->getScores();
		my $Mscore = 0;
		if( abs($ESmin) > $ESmax ){ $Mscore = $ESmin; }else{ $Mscore=$ESmax; }
		
		# pushing the max value in '@scores'
		push( @scores , $Mscore);
		
		#print("\r $i/". $#geneSets);
	}
		
	$this->{"scores"} = \ @scores  ;	
}


sub getScores($)
{
	( my $this ) = @_ ;
	
	# if it is the first call
	if($this->{"scores"} == 0)
	{
		$this->_generateRandomDistribution();
		$this->_calculateGeneSetScores();
	}
	
	return @{ $this->{"scores"} } ;
}


sub computePvalue($ $)
{
	( my $this , my $score) = @_ ;
	if( !defined($score) ){ die("argument 'score' not set !"); }
	
	# if it is the first call
	if($this->{"scores"} == 0)
	{
		$this->_generateRandomDistribution();
		$this->_calculateGeneSetScores();
	}
	
	my @scores = @{ $this->{"scores"} } ;
	
	my $counter = 0 ;
	for(my $i = 0 ; $i < @scores ; $i++)
	{
		if( abs($score) < abs(  $scores[$i]  ) ) 
		{
			$counter++;
		}
	}
	
	my $pvalue = $counter / ( $this->{"Nperm"} ) ;
	$pvalue = $pvalue / 2 ;
	
	return $pvalue ;
}


1;